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ary electromagnetic fields is given in terms of mutually uncorrelated, trans- 
versely shifted, fully coherent and polarized elementary electric-field modes. 
This representation allows one to propagate non-paraxial partially coherent 
vector fields using techniques for spatially fully coherent fields, which are nu- 
merically far more efRcient than methods for propagating correlation functions. 
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e.g., many light-emitting diodes. 
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1. Introduction 

There is a growing demand to include partial spatial co- 
herence in optical design, not only because of classical 
areas of application such as microscopy and projection 
lithography, but also because of the increasing impor- 
tance of partially coherent solid-state sources such as 
multimode lasers (including excimers) and LEDs. Meth- 
ods based on ray optics can not adequately describe all 
of the relevant issues related to coherence and polar- 
ization, which are intimately connected in electromag- 
netic coherence theory. Thus wave-optical methods are 
required, which must be capable of dealing with non- 
paraxial fields and systems that may contain also micro- 
and nanostructures in objects and interfaces. To this 
end, one needs computationally efficient physical-optics- 
based representations of spatially partially coherent elec- 
tromagnetic fields. 

Apart from some specific models that allow analytic 
solutions, propagating spatially partially coherent light 
even in free space is a formidable numerical problem in- 
volving four-dimensional integrals [l|. The dimensional- 
ity of the propagation integrals can be decreased to two 
(for planar sources) if the partially coherent field is repre- 
sented as an incoherent superposition of fully spatially co- 
herent fields. The classical way to do this is the coherent- 
mode expansion of the cross-spectral density (CSD) func- 
tion by means of Mercer's expansion which has re- 
cently been extended to electromagnetic fields [1, 0] . In 
this representation the coherent modes are uniquely de- 
fined by the CSD through a Fredholm integral equation; 
they form a complete and orthonormal set, in which the 
effective number of modes N increases as the degree of 
coherence of the field is reduced Thus, in free- 



space propagation, the original four-dimensional integral 
for CSD is replaced by N two-dimensional integrals for 
the coherent modes. In light-matter interaction analy- 
sis one solves the diffraction or scattering problem for N 
coherent fields of different functional forms [l3] • 

There is also an alternative representation of a par- 
tially coherent field in terms of uncorrelated, fully co- 
herent fields (lT| - [T5| . Here all coherent fields or 'ele- 
mentary modes' are of identical functional form but spa- 
tially (or angularly) shifted with respect to each other 
and weighted by a function determined by the CSD. Un- 
like the Mercer expansion, the shifted-elementary-mode 
representation is applicable only to a specific class of gen- 
uine CSDs 13]. However, this class contains many of the 
fields of practical significance in optical design, includ- 
ing all quasihomogeneous fields (LEDs, excimers, illumi- 
nation in microscopy and projection lithography, etc.). 
Since the elementary modes are identical, only a single 
2D integral needs to be evaluated in free-space propa- 
gation problems. In interaction problems one needs to 
scan the elementary mode across the object and perform 
a set of independent diffraction calculations for coherent 
light. Typically the elementary mode has a smooth func- 
tional form, at least compared to the higher-order modes 
in the Mercer expansion, and is therefore easy to propa- 
gate numerically. Another advantage of the shifted-niode 
model is that there is no need for numerical solution of 
the Fredholm integral equation: the elementary mode 
and the associated weight function can be determined, 
e.g., from far-zone properties of the field (Tsj . 

In this paper we generalize the scalar shifted- 
elementary-mode representation to the vectorial case, 
which is necessary to adequately model partially spatially 
coherent, partially polarized sources. There are two ma- 
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jor reasons why such a generahzation is necessary: first, 
non-paraxial fields can not be adequately described by 
a scalar model and, second, optical components in the 
system can modify the polarization state of the field. It 
turns out that the vectorial nature of the field does not 
fundamentally complicate the numerical procedure. Two 
elementary modes are needed to specify the state of po- 
larization, but also in the electromagnetic case the modes 
and their weight functions can be determined from far- 
zone properties of the field. 

We begin the discussion by briefly reviewing the scalar 
model in Sect. [2] to establish the notation and to simplify 
the interpretation of the main results. The extension to 
the electromagnetic case is outlined in Sect. [3] and the 
rigorous mathematical formulation is presented Sections 
U] and [5l The important special case of rotationally sym- 
metric and quasihomogeneous fields is discussed in Sec- 
tions ini and [71 respectively. Some numerical results are 
provided in Sect. El In Sect. |9l we apply the model to a 
simple LED geometry. Finally, issues such as the mea- 
surements required to determine the elementary modes 
and their weight functions are discussed in Sect. |TU1 

2. The scalar model 

Using the notations of Fig. [1] we may write the well- 
known relationship [l| between the cross-spectral density 
function VF(ri,r2) and the angular correlation function 
A{ki,K2) as (we omit the dependence on the angular 
frequency lu for brevity thoughout the paper) 



(a) y 



M^(ri,r2) - 



A{ki,K2) 



(2^)4 

X exp [i(k2 ■r2-kl- ri)] d^Ki d^Ks- (1) 



Here the asterisk indicates complex conjugation, Vj = 
{xj,yj,Zj) with j = 1,2 are position vectors, kj = 
{kjx,kjy,kjz) — {Kj,kjz) represent wave vectors, and 
Kj — ikjx,kj^y) are their transverse projections. 

Restricting now to the specific class of fields mentioned 
in the introduction, we assume that the angular correla- 
tion function is of the Schell-model form 

A{K,,K2)=9{AK)riKi)f{K2). (2) 

The radiant intensity of a scalar field is defined as [l| 

J(rs) = 2m:^k^ cos^ M(fccr, kcr), (3) 

where n is the refractive index of the medium, s — r/r 
with r = ||r|| is a unit direction vector, cr is its transverse 
projection, 9 is the angle between s and the z axis, and 
k = |lk|| is the wave number. Using Eq. ([2]), we have 



J(rs) = 2mr^k'^ cos^ \f{k(T)\' 



(4) 



Thus the radiant intensity of a Schell-model partially co- 
herent field defined by Eq. ([2]) is the same as that pro- 
duced by a coherent field with angular spectrum /(k). 



(b) t 




Fig. 1. Notations used for (a) position and (b) wave- 
vector coordinates. 



Let us introduce two-dimensional Fourier-transform re- 
lations 



e(r) 



1 



and 



pip) 



/(K)cxp(ik-r)d2K (5) 



5(A/t) exp (iA/« • p) d^AK, (6) 



where A/« — k,2 — k,i and p — x'k + yy is the transverse 
projection of the position vector. Inserting Eq. ([2]) into 
Eq. ([ij and using the Fourier representation of ^(Ak) 
obtained by inverting Eq. (jS]), we get the expression 

I^(ri,r2) ^ fr p{p')e{v, - p')e{v2 - p')d'p' (7) 

J J —QO 

for the CSD [l^. This result applies, in particular, at 
z = 0. Thus e(p, 0) is the coherent source-plane field with 
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angular spectrum /(k). The representation in Eq. ([7]) ex- 
presses the partially coherent field as a weighted linear 
superposition of spatially shifted but identical fully co- 
herent elementary (scalar) fields e(p, 0). 

3. The electromagnetic extension 

As is evident from Eq. ([3|) , the mathematical form of the 
elementary field mode e(p, 0) can be determined from the 
knowledge of the radiant intensity (at least apart from a 
phase factor, which can in fact be employed to model vol- 
ume sources jll|). In general, the weight function p{p') 
can be determined from far-field coherence measurements 
using Eqs. ^ and ©J but often there are simpler ways to 
at least approximate it fl3| . This is the case, in particu- 
lar, if the field is quasihomogeneous, i.e., if the coherence 
area in the source plane is much smaller than source area; 
in this case the weight function comes out of the integral 
in Eq. ([7]). When e(p, 0) and p{p') are known, coherent 
propagation techniques for e(p, 0) and a linear superpo- 
sition according to Eq. ([7]) suffice to propagate the entire 
spatially partially coherent field. 

The extension of the scalar shifted-mode model to the 
electromagnetic case is not trivial. One might be tempted 
to use some superpositions of, e.g., locally linearly polar- 
ized modes of the scalar functional form to model par- 
tially polarized or unpolarized sources. However, such 
constructions seem hard to justify mathematically. The 
approach taken here is based on far-field information: in 
the far zone the field behaves as an outgoing spherical 
wave, and therefore it has a well-defined local polariza- 
tion state. We employ this fact to separate the angular 
correlation tensor, which is the electromagnetic exten- 
sion of the function A(ki, K2), into two orthogonal parts 
[see Eq. ([50]) in sect.|3]. These represent the electromag- 
netic elementary modes (or polarization modes) in the far 
field. Then the source-plane modes can be determined by 
Fourier-transform techniques in analogy with the scalar 
case. 

As a result of the construction process (to be described 
mathematically in the following sections), we obtain a 
vectorial shifted-mode expansion for the spatially par- 
tially coherent field everywhere in space [see Eqs. (15^ . 
([33| . and (|40l) in Sect. [5]. Instead of propagating one 
fully coherent mode as in the scalar case, we now need 
to propagate two well-defined fully coherent vectorial 
field modes, and to form the generalized shifted-mode 
superposition, to govern the propagation of the electro- 
magnetic spatially partially coherent field. Thus the 
increase in computational complexity, compared to the 
scalar case, is essentially fourfold. 

4. Field representation in the far-zone 

A statistically stationary random electromagnetic field 
in the space-frequency domain is described by a CSD 
matrix W(ri,r2), which may be expressed as (I. ITg} 



Here T indicates the transpose, the brackets denote en- 
semble averaging, and the electric-field realizations E(r) 
are understood as appropriate random linear superpo- 
sitions of the eigenfunctions of the Fredholm integral 
equation satisfied by W(ri,r2) [1, |3l- The position- 
dependent spectral density of the field can be written, 
in analogy with scalar theory of partial coherence Q , as 
S{y) — tr W(r, r), where tr stands for trace. 

The relation between the field at the (secondary) 
source plane z = and the far-field can be found, for 
example, using the angular spectrum representation of 
the CSD matrix 17]. We thus have, at any plane z > 0, 



W(ri,r2) 



(27r)4 

X exp [i(k2 • r2 - kj • n)] d^Ki d^K2, (9) 



where the angular correlation matrix (ACM) 

A{Kun2)= JJJJ W(Pi,P2,0) 

X exp [i(Ki • pi - K2 • P2)] d^pi (fp2 (10) 

describes the correlations between vectorial plane-wave 
components. In the far zone [l7l | 



W°°(riSi, r2S2) = (27rfc)^ cos 01 cos92A{kcTi, kcT2) 
exp [ifc(r2 — ri)] 



rir2 

and the Poynting vector takes the form [l3] 

P°°(rs) = !J &rW°°(rs,rs) 
2 V Mo 



(11) 



= § cos 



2r7,7r2fc2 



l — trAik(T,k(T), (12) 
Mo 



where n is the refractive index of the material, and eg 
and fiQ are the vacuum permittivity and permeability, 
respectively. Furthermore, the radiant intensity is 



J(rs) = lim p||P°°(rs)||] 

r— >oo 

= 2mT^k'^ cos^ ej—tr A{k(T, fccr) 



(13) 



Owing to Eqs. and pUj) . also ACM has a represen- 
tation as a correlation matrix: 



A{k^,K2) - (A*(/«i)AT(k2)), 



(14) 



where the components of A(k) represent, component- 
wise, the angular spectra 1] of the electric-field realiza- 
tions. It follows from Eq. ([Ti|) that A{n, k) is Hermitian, 
satisfying 



W(ri,r2) = (E*(ri)ET(r2)). 



(8) 



A\ki,K2) = A{k2,Ki), 



(15) 
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where the dagger denotes the adjoint matrix. It is also 
non-negative definite in the sense that 



SL'{Ki)A{Ki,K2)a{K2)d^Kid'^K2 > 0, (16) 



where sl{k) is an arbitrary, sufficiently well-behaved vec- 
tor function of the same size as the electric-field vector. 

The following argument is essential for the conclusions 
of this paper: the field in the far zone is well known 
to be a modulated outgoing spherical wave. Thus, in 
spherical polar coordinates, the angular-spectrum vector 
is two-dimensional, i.e., s • A(fccr) = for every s. Hence 
A{ki, K2) is expressible as a 2 x 2 matrix in these coor- 
dinates. 

Since A{ki,K2) is a Hermitian, non-negative definite 
2x2 matrix, it has two non-negative real- valued eigen- 
functions. In view of Eqs. (fTS)) and ([16]), we have at 
Ki = k,2 = K, the decomposition 



With the aid of Eqs. ^7} and ([18]), we then obtain 



P(rs) 



h{k(T) - hika) 



hika) + hika) 



(22) 



Thus P(rs) = 1 for each individual polarization mode, 
i.e., if either Ii{k) = or /2(/«) — 0. It is worth stressing 
that the polarization decomposition presented above is 
analogous with the well-known decomposition of a par- 
tially polarized plane wave into two polarization modes: 
in particular, Eq. (|22]) is analogous with Eq. (6.3-31) in 
Ref. 111. In our case, however, the direction-dependent 
eigenvalues are those of a partially polarized and par- 
tially coherent field in the far zone. 

To be able to include partial coherence (in addition to 
partial polarization) in the analysis, we now assume that 
also the two-point ACM can be expressed in the form 



A{ki,K2) = Ai{ki,K2) +A2{Ki,K2) 



(23) 



^(k, /«) — /j(/t)Fj (k.)Fj (k.), (17) where Aj{Ki,K2) have diagonal values Aj{K,K) = 



where /j(/«) are the eigenvalues and Ff{K) are the eigen- 
vectors of A{k,k). The eigenvectors may be assumed 
orthonormal, i.e., 



(18) 



In other words, the matrix A{k, k) can be diagonalized. 
Evaluation of the trace in Eq. (|T^ yields the explicit form 
for the radiant intensity: 



J(rs) = Jo cos^ [/i(fccr) + hika-)] 



(19) 



where Jq = 2mr'^ k'^ eo / . Thus it is proportional to 
the sum of the eigenvalues of the ACM at Ki = k,2- 

The decomposition in Eq. ([TT]) is, in fact, valid regard- 
less of the chosen coordinate system; the eigenvalues re- 
main invariant in all unitary transformations, including 
simple coordinate transformations between, e.g., Carte- 
sian and spherical polar coordinate systems. As explicitly 
expressed in Eqs. (jTT]) . the polarization decomposition of 
ACM is generally direction-dependent: the eigenvalues 
and/or the eigenvectors depend on s. 

The physical meaning of Eq. ([TT]) is clear: In each di- 
rection in the far-zone, we may decompose the single- 
point ACM into two mutually uncorrelated, orthogonal 
polarization components. Moreover, owing to the factor- 
ized form /j(K:)F*(K:)Fj(«:) of these components, they 
both represent fully polarized fields. This can be verified 
from the well-known formula for the (space-frequency 
domain) degree of polarization: 



/j(K)F*(K)Fj(K). In other words, A{ki, K2) is assumed 
to be expressible as a sum of two mutually uncorre- 
lated, but fully coherent and polarized modes even if 
Ki 7^ K2- We stress that, while Eq. ([^^ does not fol- 
low from Eq. (|17p . it typically holds. We do not dwell 
into a detailed discussion of this point here, but note that 
it is a nontrivial task to find counterexamples. 

We proceed to investigate the angular correlation prop- 
erties of the class of fields described by Eq. ([^5]) . Without 
loss of generality, we may employ spherical polar coordi- 
nates in the far field, in which case Aj{Ki, K2) is a 2 x 2 
matrix. Let us denote by U{k) a unitary matrix whose 
columns are the eigenvectors Fj(/t) of A{k,k) and by 
^{k, k) a diagonal matrix with elements equal to the 
eigenvalues Ij{K) of A{k,k). Then, equivalently with 
Eq. (|17|). we may write 



A{k, k) = U*{k)T>{k, k)U^{k) 
and consequently 

U^{k)A{k, k)U*{k) = T>{k., k). 



F(r) = 1 - 



4dctW(r,r) 



1/2 



[tr W(r,r)]2 
In the far zone we have, using Eq. (|lip. 

4 det ^(fccr, fccr) ^ 
[trA{k(T,kcT)]^ 



P(rs) = <l- 



(20) 



(21) 



(24) 



(25) 



In view of Eq. ([23]) and the associated discussion, 
Aj{iii, K2) represents the angular correlation between 
two plane- wave components whose polarization states are 
described by deterministic vectors Fj(/«i) and Fj(/«2)- 
Thus we may write, in analogy with Eq. (|17p . 

2 

A{ki,K2)=J2Gj{ki,K2)F*{ki)FJ{k2), (26) 
i=i 

where Gj{Ki, K2) are scalar (angular correlation) func- 
tions. It is seen by direct calculation that A{k,i, K2) has 
a representation similar to Eq. ([^5]) . 

U^iKi)A{Kl,K2)U*iK2) = X»(ki, K2), (27) 
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where X>(ki,K2) is a diagonal matrix with elements 
Gj{Ki,K2)- Thus Fj(/c) and Gj{Ki,K.2) can be deter- 
mined by singular value decomposition of K2). 

It can be shown, e.g., using Eq. (|16l) that the singular 
values Gj{K.i,K2) satisfy 



|G,(/^l,/^2)|'</,(Kl)/,(/^2). 



(28) 



As a result, we may define the normalized angular corre- 
lation functions 



Gj{Ki,K,2) 
[I,{K,)I,{K2)f' 



(29) 



satisfying the inequalities < \gj{Ki,K2)\ < 1. If the 
correlations in the far zone are of (generalized) Schell- 
model form, i.e., gj{Ki,K2) — gj{A.K), it follows from 
Eqs. (HH) and ^ that 

2 

A{kuK2) - ^g,(A/^)f;(Ki)f7(/^2) (30) 

with 

f,(Ac) = [/,(K)]V2F,-(/.). (31) 

This is the electromagnetic extension of Eq. 
5. Elementary electric-field modes 



It follows from Eq. ([23]) and the linearity of Eq. ([TOj) that 
the CSD has the decomposition 

W(ri,r2) = Wi(ri,r2)+W2(ri,r2). (32) 

Let us define the inverse Fourier transforms ej(r) and 
Pj{p) of the functions fj(/«) and gj{/S.K) in analogy with 
Eqs. ([SJ and (jH]). With a procedure similar to that used 
in derivation of Eq. ([7]), we can express the two terms in 
Eq. ([32]) in the form 



W,(ri,r2) 



p,(p')e*(ri-p')eJ(r2-p')dV' 



(33) 



in the half-space z > 0. In particular, Eq. ((33|) is valid 
at the source plane 2 = 0, where ej{v) = e_,(p, 0). 
This result is the electromagnetic extension of the scalar 
elementary- mode decomposition in Eq. ([7]). 

Let us next examine some general properties of 
Eq. (1221). Equations (HH]) and §Bl together with 

fji^) =11 ej (r) exp (-ik • r) d^p (34) 

lead to 



el(p-p',0)e2(p,0)exp(-i/«-p') d^pd^' = 0. 

} 

(35) 



Because of the exponential factor in the integral, ei (p, 0) 
and e2 (p, 0) are generally not orthogonal in a point- 
wise sense. For many paraxial fields, though, the 
functions fi(K) and f2(/«) are globally orthogonal, i.e., 
f|(Ki)f2(/«2) = for all Ki and K2. In this special case 
the pointwise orthogonality of ei(p, 0) and e2(p,0) fol- 
lows from Eqs. ([TH]), and This property holds 
also for significant classes of non-paraxial fields, as will 
be demonstrated in Sections |6] and [H However, it is not 
essential for practical implementation of the propagation 
algorithm developed in this paper. 

Since the functions e^ (r) represent fully coherent fields, 
they obey the Helmholtz equation 



(36) 



Furthermore, the CSD matrix obeys the divergence equa- 
tion (4|] 



VTw(ri,r2) = 0, 



(37) 



where the subscript 1 denotes differentiation with respect 
to ri. Together with Eqs. (15^ and ([25)) . this implies that 



V-e,(r)==0. 



(38) 



In view of Eqs. and ([25)) . the functions fj(r) behave 
exactly as the electric field vector in the space-frequency 
domain. Hence the matrix-functions 



WKri,r2)=e*(ri)eT(r2) 



(39) 



may be called the elementary electric- field modes of the 
field. Using these modes, we can write Eqs. ([5^ and ([25)) 
in a more compact form 

W(ri,r2) = ^ // p,(p')W,^(ri-p',r2-p')dV 



and express the spectral density of the field as 



(40) 



5W = E 



p,{p')S^Av~p')A^p', (41) 



where S'l(r) = trW^(r,r) are the spectral densities of 
the two elementary electric-field modes. 

In analogy with scalar theory [l3] , the field in Eq. (|in)) 
is understood to consist of a weighted continuum of iden- 
tical, laterally shifted elementary modes. The main dif- 
ference between the scalar and electromagnetic descrip- 
tions is that the electromagnetic field consists of two (un- 
correlated) sets of elementary modes, whose polarization 
states in the far-zone are orthogonal and whose radiant- 
intensity distributions are in general different. Owing to 
the factorized form of Wj(ri,r2) in Eq. ([2^ . each ele- 
mentary mode is completely coherent [18j in the sense of 
the space-frequency analog of the degree of coherence for 
electromagnetic fields put forward in Ref. [l^ : 



/^(ri,r2) = {tr[/x(ri,r2)p(r2,ri)]} 



1/2 



(42) 
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Here 



6. Rotationally symmetric fields 



M(ri,r2) = 



>V(ri,r2) 



1/2 



(43) 



is the normalized CSD matrix. In general, all 3 x 3 ma- 
trix elements of W°(ri,r2) are non-zero and spatially 
varying. 

Analogously with Eqs. (H^ and ([^S]). we can define the 
degree of angular coherence 



a{Ki,K2) = {tr [a{Ki,K2)a{K2, Ki)]} 



1/2 



where 



a(Ki, ^2) 



[tv A{Ki, Ki)tT A{k2, K2)] 



1/2 



(44) 



(45) 



is normalized angular correlation matrix. 

To conclude this section we provide a convenient, fully 
general series representation for the electric-field modes. 
In the far zone, where the field is transverse to the local 
propagation direction and thus fluctuates locally in the 
plane defined by the spherical polar unit vectors 6 and 
i/j, it is natural to express the basis vectors in the form 



(46) 



Using spherical polar coordinates {0,ip) for the spatial 
frequencies and circular cylindrical coordinates (p, 0, z) 
for the position vector (see Fig. [T]), we then have (after 
somewhat lengthy calculations outlined in Appendix |A|) 



— i™ exp(iTO0) / exp (ikz cos 9) 

/n 

rvi — rv^. 

d 



m— — 00 

m 



dp 

Tfi d 
— COS 9 fj^e^„i{9) - ifj.4,,m{9)-^ 



zksm^9f,^eM0) } Jm{kp sin 9) d9, (47) 



where 



2-R 



hi,m{S) = TT / hd^' ^) (-i™^) dV' (48) 

are the azimuthal Fourier coefficients of fj,^{9,ip) and 
^ stands for either 9 or ip. Note that the upper limit 
of the integral in Eq. (|47|) is set to 7r/2. Thus the el- 
ementary electric-field modes are taken to contain only 
propagating waves, i.e., information that can be gathered 
from far-zone measurements. As a result, the propaga- 
tion method considered here is not suitable for modeling 
near-field phenomena (fields at distances of the order of 
one wavelength from the plane z = 0). 



Let us assume that the CSD at 2 = is rotationally 
symmetric about the z axis. Then also the ACM is ro- 
tationally symmetric and the polarization basis vectors 
fi(/«) and f2('«) are rotationally invariant, i.e., their 9 and 
■0 components fjfi{9,ip) and fj^^{9,ip) are independent 
on the azimuthal angle ^. Hence only the zeroth-order 
Fourier coefficients fj,^fl{9) in Eq. P5|) are non-zero and 
it follows from Eqs. ^\ and (|A14p that 

e, (r) = — / sin 9 exp (ikz cos 9) 
Jo 

X <^ iJi(A:psine') cos 6'/j_e_o(6')p + /j,^,o(^)0 



sin9Jo{kpsin9)fj^eA^)^ \ dfi* 



(49) 



Thus, as expected, the elementary electric-field modes 
are also rotationally symmetric about the z axis. In par- 
ticular, if the basis vectors are parallel to 9 and 

1,2 /•7r/2 

ei(r) = — / /..M^sin^ 

Jo 

X [i/5Ji(fcpsin6') cos 6* — z sin0Jo(fc/3sin6')] 
X exp (ifcz cos (?) d9, (50a) 

e2(r) = 0— / f j^^fl{9) sin 9 J i{kp sin 9) 

X exp {ikz cos 9) d9. (50b) 

The elementary modes ei(r) and e2(r) are now radi- 
ally and azimuthally polarized fields, respectively. Hence 
they are pointwise orthogonal, regardless of whether the 
field is paraxial or not. Since these modes have no az- 
imuthal phase variation, they possess no phase singular- 
ity (vortex) at p = 0. Thus, even if the axial field van- 
ished z = (as turns out to be often the case for some of 
the field components), such a zero does not propagate. 

7. Quasi-homogeneous sources 

Assume next that the variations of the spectral density 
S{p, 0) at the source plane are slow compared to the vari- 
ations of the degree of coherence ^(p^, 0, P21 0)- More- 
over, let the correlations at z = be of the Schell-model 
form, i.e., depend only on Ap = P2 ~ Pi- Such a pla- 
nar source is said to be quasi- homogeneous We may 
then approximate S{pi,0) ~ S{p2,0) ~ S{p,0), where 
p ~ (pi + P2)/2, and write 



W(Pi,P2,0)«5(p,0)MAp,0). 



(51) 



Inserting Eq. ([51]) into Eq. ([TU]) and defining k — (ki 
K2)/2 yields 



(52) 
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S{p, 0) exp {-iAk ■ p) d^p, (53) 



/x(Ap,0)cxp(-iK • Ap)d^Ap. (54) 



where 
and 
Since 

tr A{k,k) = S{0)trfi{K) 

= 5(0) ff tr/x(Ap,0)exp(-iK • Ap)d2Ap, 

(55) 

it follows from Eq. ()13p that the radiant intensity pro- 
duced by a quasihomogeneous electromagnetic source de- 
pends only on the correlation properties of the source 
field, in complete analogy with the scalar case [l]. 

Inserting from Eq. ((5T|) into Eq. pS)) and recalling that 
/i is a wide function compared to S, we have 



q:(ki, K2) 



S{Ak) /2(k) 

S{0) Ko)' 



(56) 



Using Eq. (|44| and noting that is real because S is 
non-negative, we obtain 



SiAK) {t"" [m('^)]'} 



1/2 



5(0) 



tr /t(/«) 



(57) 



This expression can be cast into a more transparent 
form using the far-zone degree of polarization defined in 
Eq. (j2ip . which can be written equivalently in the form 



P{k) = 



2tr [/i(«)] 



1/2 



[tTfl{K)Y 

We then have, from Eq. ([ST)) . 



- 1 



a(Ki,/«2) = 



5(0) 



1/2 



(58) 



(59) 



The first fraction in this expression is equal to the scalar 
degree of angular coherence. The second fraction, how- 
ever, is a polarization-dependent modulating term that 
depends on the source-plane correlations. 

The assumption that the source is quasihomogeneous 
simplifies decisively the elementary-mode decomposition 
of the field in the scalar case [l^ , and the same is true in 
the electromagnetic case. It follows from Eqs. p7)) . (I5T]) . 
and mi that 



E// e*(pi-p',0)ej(p2-p',0)d 
A-{n, k) exp(z/^ • Ap) d^/^ 



(2^)2 
5(0)Ai(Ap,0), 



(60) 



where, in the last step, we have used Eq. (1521) . Comparing 
Eqs. dSlI) and ^ yields 

5(p,0) 



>V(Pi,P2,0) 



5(0) 



2 

E 



e*(pi-p',0)eKp2-p',0)dV. (61) 



Thus the weight function no longer appears inside the in- 
tegral and, irrespective of the spatial distribution of the 
spectral density at the source plane, the field character- 
istics can be determined the propagating a convolution 
integral involving the elementary modes only. 

8. Sources with cosine-power radiant intensity 

Let us consider the rotationally symmetric case with radi- 
ally and azimuthally polarized basis vectors Fi (0, -0) = 0, 
F2(^, V") = '01 a-nd corresponding eigenvalues Ii{9,ip) — 
A\ cos°"^ 6*, /i(6i, V') = Al cos''"^ 9, where Ai and A2 are 
arbitrary (real) functions of frequency. Then, in view 
of Eq. (jl9p . the radiant intensity is a superposition of 
two cos" 9 type contributions, one radially and the other 
azimuthally polarized: 

Jo [Al cos" 9 + Al cos*" 9] . 



j(0,V) 



(62) 



The elementary electric-field modes in the far zone are, 
according to Eq. (|5T|) . 



fi(6',?A) = OAi cos°'''^-^9, 
f2(6',V) =i/;A2 cos^'^-^9. 



(63a) 
(63b) 

Using Eqs. and we see that the non- vanishing 
Fourier coefficients are fi^e.Q — Ai cos"^/^"^ 9 and /2.v;.o — 
A2 cos*"/^"^ 9. Inserting these into Eqs. ([50]) and apply- 
ing (IB4bl) derived in Appendix [B] we obtain the source- 
plane elementary field modes in the form 



ei(p,0) = Al 



X 1F2 

- zT 



^ 2 



1F2 



1 

2 + 
a 

3.1- 

2' ' 2 



e2(p,0) 



^ ^16V7F V4 



X 1F2 



^2^ 
2' '2 



b 

4'" 



a 
4 

4 

3 a 
~ ^ 4' 



(64a) 



(64b) 



where F is the Gamma function and 1F2 is the regular- 
ized hypergeometric function (see Appendix [B]) . 

Figures [SHU illustrate the radial dependence of the 
elementary-field components and the function 



w{p,0) ^ [\\ei{p,0)r + \\e2{p,0)\\T 
for different values of a = 6 = n, with Ai — A2 



(65) 

-ifc-2. 
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Fig. 2. Relative amplitudes of the radial (solid line), 
azimuthal (dashed line), and longitudinal (dotted line) 
components of the elementary fields as a function of the 
normalized radial coordinate kp, as well as the function 
w{kp) (thick solid line) for ti = 1, which corresponds to 
a Lambertian source. 
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Fig. 5. Same as Fig. |4l but with n = 20. The radiation 
pattern is increasingly directional and could be produced 
approximately by a LED with an integrated collimating 
lens. 

9. Illustration: LED model 
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Fig. 3. Same as Fig. [31 but with n = 2, which corresponds 
to an incoherent source in scalar theory. 
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Fig. 4. Same as Fig. [21 but with n — 5. Thus the source 
has a somewhat directional radiation pattern. 



Let us consider a simple model for a rotationally symmet- 
ric surface-emitting LED illustrated in Fig. [6^, where the 
primary light-emitting region is planar (such as a quan- 
tum well) and buried inside a semiconductor material of 
refractive index rig. Each primary source point is as- 
sumed to radiate (independently) a spherical wave. If we 
denote the propagation angle inside the semiconductor 
material by 9' , the radiant intensity may be expressed as 
a sum of radially and azimuthally polarized contributions 



J. 



J = 1,2: 



7W/ 



r) = cos^ 9' I. 



(0/ 



r) = j«, cos^ 



A 



(0/ 



(66) 

where A^'^ {9') is the complex amplitude of the plane wave 
in direction 9'. If the distance between the pn plane 
and the semiconductor-air interface is large compared to 
A, the local plane wave approximation is valid at the 
semiconductor-air interface. Thus the output radiant in- 
tensity takes the form 



Jj{9) = Jo J cos^ 91^ {9) = Jo J cos^ 9 \Aj{9)y 



(67) 



Here Joj = Jqj/us, the angles 9 and 9' are related by 
SneU's law sin 6* = nsSm9', and Aj{9) = tj{9,9')Aj{9'), 
where tj{9, 9') are given by Fresnel's equations 



2ns cos 9 



cos 9' + n.s cos 9 ' 



t2 



2 cos 6* 



Us COS 9' + cos 9 



(68) 



(69) 



for radial (TM) and azimuthal (TE) polarizations. 

Because of the large refractive index Ug of a semicon- 
ductor, only a narrow cone of plane waves with incident 



angles 9' in the range < 0' < arcsin (1/ns). We may 
thus assume that the radial and azimuthal contributions 
to the radiant intensity of the primary source are equal, 



(a) 



I.e. 



_ 7-(i) 



and hence we may denote Jq = 2 J( 



and Jo = 2 Jo J . Then the degree of polarization given by 
Eq. ^ is 



p{e) = 



(70) 



and the radiant intensity transforms at the interface ac- 
cording to 



J{e) _ 1 cos^6> 
JW(6») ~ 2ns cos2 9' 



\h{9,9'f + M9,9'f 



(71) 



If we assume that the radiation pattern produced by the 
primary source is Lambertian, with 



jf{9')^^-f^>cos9', 



1 



(72) 



the radial and azimuthal contributions to the radiant in- 
tensity of the secondary source become 



(73) 



These contributions Ji{9) and J2{9) are shown by the 
dotted and dashed lines, respectively, in Fig. ^p, where 
we have taken rig = 3.5. The curves fit well the cos" 9: 
we obtain n ~ 3.4 = h for the azimuthally polarized 
contribution, n w 2.4 = a for the radially polarized con- 
tribution, and n ~ 2.9 for an equally weighted sum of 
the two contributions. Thus the elementary electric-field 
modes given by Eq. (j64p provide good approximations of 
the modes of the structure in Fig. [6^. 

The degree of polarization, also plotted in Fig. ^p, 
increases from a zero on-axis value (unpolarized radi- 
ation in the paraxial domain) to P{9) ~ 0.85 when 
9 — >■ 7r/2, indicating partially polarized radiation in the 
non-paraxial domain. 

10. Final remarks 
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Fig. 6. (a) A generic geometry of a broad-area surface- 
emitting LED: pn is the active emitting area and n is 
the refractive index of the semiconductor material. The 
solid and dashed curves illustrate the azimuthally and ra- 
dially polarized contributions to the radiant intensity dis- 
tribution J {9). (b) Geometrical-optics predictions of the 
azimuthally (dashed curve) and radially (dotted curve) 
polarized contributions to the radiant intensity, their av- 
erage (solid curve), and the degree of polarization P{9) 
in the far zone (thick solid curve). 



the field in the far zone. Thus only single-point mea- 
surements across the radiation pattern are needed. De- 
termination of the weight functions requires, in general, 
two-point correlation measurements in the far zone, but 
this is avoided if the field is known to be quasihomoge- 
neous. 



The electromagnetic elementary-mode decomposition 
presented in this paper should prove useful in in opti- 
cal system modeling by field tracing methods. To this 
end, it is necessary to determine elementary field modes 
and the weight functions of the source. The example 
presented above illustrates the possibility of doing this if 
there is sufficient a priori information about the structure 
of the source. If, however, the source properties are not 
known, it is necessary to determine the modal decompo- 
sition experimentally. As in the scalar case [Tsj , this can 
in principle be accomplished by far-field measurements. 
In general, the polarization basis vectors fj and the eigen- 
values Ij can be determined from the polarization matrix 
A.{k,k). This matrix can be determined by measuring, 
e.g., the angular dependence of the Stokes parameters of 
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Appendix A: Derivation of Eq. (47) 

In this Appendix we present some details of the deriva- 
tion of Eq. (j47|) . which is the general representation of the 
elementary electric-field modes. Denoting the circular 
cylindrical coordinates by (p, z) in the position-vector 
space and by (k, ip, z) in the wave vector space, and the 
spherical polar coordinates by (fc, 0, ij)) in the wave vector 
space, we have the following relations between the unit 
vectors of these systems and the Cartesian coordinates 
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(see Fig.[T]): 



= R(0) 



R(^) 



ll{0) 



where 



cos f sin ^ 
-sin^ cos^ 



(Al) 



(A2) 



(A3) 



(A4) 



is the rotation matrix and ^ may stand for cf), tp, or 9. 
Using Eq. (|A3I) and recognizing that fcz = kcosO and 
K = fcsin^, we can write Eq. (1461) in circular cyhndrical 
coordinates: 

ijin) cos6fj^g{K)k + fj^.^{K)'ip - sin 61/^,0 (k)z. (A5) 
Using Eqs. (|Aip and (jA2l) one can estabhsh the relation 





1 


■ 1 1" 




~ 2 


— i i 



3''^ [p + i-0^ 



(A6) 



where /? 
f,(0,^) 



7/;, and express fj(/«) in the form 



- [cos ef^,e{e, tA) - i/j, ^ (0, V)] c''' (p + it^ 

i [cos 0/,, e((?, tA) + i/,,^ (0, ^)] e-'^ (p - it^ 



- sin6'/j,e(6'»z 



(A7) 



Recognizing further that k-r — n ■ p + kzZ and that 
K - p — kp sin 9 cos /3, we can expand the plane- wave term 
in Eq. ([6]) in a series form 



exp (ik • r) = cxp (ifcz cos 6') i™ Jm(fcpsin0)e 



im/3 



(A8) 

where denotes the Bessel function of the first kind 
and order m. Combining Eqs. (??) and (|A7|) we then 
have 

= i" / gj,m(r)exp(ifc0cos6') 

m=-oo 

X J„i(fcpsin6')fcsin6'd6', (A9) 



where 



(AlO) 



Performing the integration with respect to with the aid 
of the definition of the azimuthal Fourier coefficients in 
Eq. (jlH]), we obtain 

gj,m(r) i (^p + i0^ [cos6'/j,9,,„+i(6') - i/,-^V':"i+i(^)] 

xe'(™+i)'^J„(fcpsin6l)fcsin6' 

+ ^ (^p - i<^^ [cos9fjfi^m-i{9) + i/j,^,„_i(6')] 

xe'(™-i)"*J„(fcpsin6')fcsin6l 
- ze™'^J™(fcpsin6')fcsin2 6l. (All) 



Collecting terms in the summation of Eq. (IA9|) by re- 
placements 

e'("±i)'^J„,(fcpsin6') ^ i=Fie""'^J™^i(fcpsin6') (A12) 
and employing the Bessel-function identities [l^ 

2m 

J™_i(x) + Jto+i(x) = — Jm(a;) (A13) 

X 



Jm~l{x) - Jm+l{x) ^ 2 — Jm{x) (A14) 



and 



we finally arrive at Eq. (|47)) . 

Appendix B: Derivation of an integral formula 

Let us first recall the series-representation of Bessel func- 
tions lii: 



^^(7) = E 



j!(m + j)! 



(i) 



On the other hand, we have the relation [2l| 



sm^ 9 cos^ 9d9 = -B , 

2 V 2 2 



where 3?(p) > -1, ^{q) > -1, 



r(7 + 



(Bl) 



(B2) 



(B3) 



is the beta function, and r(7) is the Gamma function [20] . 
Combining Eqs. (|B1[) and (IB2p . and interchanging the 
order of summation and integration yields 



-it/2 



sinP 9 cos'' 9 J,n{kp sin 9) d9 



{kpT fq + l 

2m+l 



^[j + Up + m + 1)] 



Y 

^ j\r{m + j + i)r[j + 1 + i(p + TO + q)] 



2^2 \ J 



-k^p 



(B4a) 
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where we have employed the identity s! = T(s + 1) for 
integer s. Equation (|B4aP can also be expressed in the 
form 

7r/2 



smPQcos'^ 9 Jmikp sine) d0 



2™+i V 2 

^ /p + m + 1 
X 1F2 



m + 1, 1 



p + q + m —k^p^ 



(B4b) 



where ^ (ai, 02, . . . , a^; &i, 62, . . . , 6<j; 7) denotes the 
regularized hypergeometric function, defined by 



pFq {ai,a2, . . . ,ap;bi,b2, . . . ,bg;-f) 
_ pFq (ai, 02, ■ ■ ■ , Qp; bi,b2, ■ . . , 7) 

r(6i)r(62).-.r(M 



(B5) 



pFq (oi, 02, . . . , flp; 61, 62, . . . , bq-, 7) is the generalized hy- 
pergeometric function [22] 



pFq (ai, 02, . . . , Op; 61, &2, • ■ • , &g; 7) 
(ai)i(a2)j-.-(ap)j7^ 



E 



and 
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